Сcылка на практикум

Целью данного занятия является ознакомиться с возможностями докинга низкомолекулярного лиганда в структуру белка.

In [2]:
from xmlrpclib import ServerProxy
from IPython.display import Image
import os, sys

# pymol launching
import __main__

__main__.pymol_argv = [ 'pymol', '-x' ]

### Если вывод в графическое окно тормозит или не нужен, то:
##__main__.pymol_argv = [ 'pymol', '-cp' ]


import pymol
 
pymol.finish_launching()
 
from pymol import cmd

from IPython.display import Image

def Img(name):
    cmd.ray()
    cmd.png(name)
In [3]:
%%bash
export PYTHONPATH=${PYTHONPATH}:/home/preps/golovin/.local_numpy/lib/python2.7/site-packages:/home/preps/golovin/.local_numpy/lib/python2.7/dist-packages
export PATH=${PATH}:/home/preps/golovin/progs/bin
Докинг одного из мономеров сахара (NAG) из прошлого занятия

Запишем SMILES-нотацию N-ACETYL-D-GLUCOSAMINE (NAG) в файл nag.smi:

In [2]:
%%bash
echo "CC(=O)N[C@H]1[C@H](O)O[C@H](CO)[C@@H](O)[C@@H]1O" > nag.smi

Построим 3D структуру этого сахара в pdb формате:

In [4]:
%%bash
obgen nag.smi > nag.mol
babel -imol nag.mol -omop nag.mop -xk "PM6"
MOPAC2009.exe nag.mop
babel -imopout nag.out -opdb nag.pdb

A T O M   T Y P E S

IDX	TYPE
1	1
2	3
3	7
4	10
5	1
6	5
7	1
8	5
9	6
10	6
11	1
12	5
13	1
14	6
15	1
16	5
17	6
18	1
19	5
20	6
21	5
22	5
23	5
24	28
25	21
26	5
27	5
28	21
29	21
30	21

F O R M A L   C H A R G E S

IDX	CHARGE
1	0.000000
2	0.000000
3	0.000000
4	0.000000
5	0.000000
6	0.000000
7	0.000000
8	0.000000
9	0.000000
10	0.000000
11	0.000000
12	0.000000
13	0.000000
14	0.000000
15	0.000000
16	0.000000
17	0.000000
18	0.000000
19	0.000000
20	0.000000
21	0.000000
22	0.000000
23	0.000000
24	0.000000
25	0.000000
26	0.000000
27	0.000000
28	0.000000
29	0.000000
30	0.000000

P A R T I A L   C H A R G E S

IDX	CHARGE
1	0.061000
2	0.569000
3	-0.570000
4	-0.730100
5	0.300100
6	0.000000
7	0.560000
8	0.000000
9	-0.680000
10	-0.560000
11	0.280000
12	0.000000
13	0.280000
14	-0.680000
15	0.280000
16	0.000000
17	-0.680000
18	0.280000
19	0.000000
20	-0.680000
21	0.000000
22	0.000000
23	0.000000
24	0.370000
25	0.400000
26	0.000000
27	0.000000
28	0.400000
29	0.400000
30	0.400000

S E T T I N G   U P   C A L C U L A T I O N S

SETTING UP BOND CALCULATIONS...
SETTING UP ANGLE & STRETCH-BEND CALCULATIONS...
SETTING UP TORSION CALCULATIONS...
SETTING UP OOP CALCULATIONS...
SETTING UP VAN DER WAALS CALCULATIONS...
SETTING UP ELECTROSTATIC CALCULATIONS...

S T E E P E S T   D E S C E N T

STEPS = 500

STEP n       E(n)         E(n-1)    
------------------------------------
    0     577.792      ----
   10    111.36434    113.84533
   20    97.45839    98.39372
   30    90.91881    91.35152
   40    87.42757    87.58316
   50    85.53617    85.67561
   60    84.30098    84.40171
   70    83.39829    83.47296
   80    82.68933    82.75073
   90    82.10277    82.15595
  100    80.31998    81.57788
  110    79.57716    79.60921
  120    78.16647    78.17204
  130    75.43072    75.87105
  140    75.11357    75.16446
  150    74.30726    74.33649
  160    74.11937    74.13423
  170    73.93195    73.94572
  180    73.75789    73.76667
  190    73.57861    73.61573
  200    73.43041    73.44388
  210    73.25618    73.29716
  220    73.10436    73.11853
  230    72.92128    72.96391
  240    72.75808    72.77292
  250    72.55725    72.59991
  260    72.38003    72.39134
  270    72.18184    72.19937
  280    71.98090    71.99282
  290    71.69521    71.71213
  300    71.43131    71.44613
  310    70.43268    70.43749
  320    70.05117    70.13624
  330    69.77147    69.78139
  340    67.78301    67.82620
  350    66.73414    66.77115
  360    66.45726    66.47988
  370    65.56069    65.58031
  380    65.32462    65.34350
  390    64.60085    64.61862
  400    64.39423    64.42835
  410    64.25034    64.25922
  420    64.10583    64.11797
  430    63.97292    63.98536
  440    63.82470    63.83571
  450    63.70571    63.71660
  460    63.49249    63.58417
  470    62.94184    62.94650
  480    62.59190    62.76601
  490    62.46526    62.47069
  500    62.24316    62.24885

W E I G H T E D   R O T O R   S E A R C H

  NUMBER OF ROTATABLE BONDS: 3
  NUMBER OF POSSIBLE ROTAMERS: 72
  INITIAL WEIGHTING OF ROTAMERS...

  GENERATED 250 CONFORMERS

CONFORMER     ENERGY
--------------------
     1        59.921
     2        58.982
     3        62.230
     4        58.893
     5        58.753
     6        59.043
     7        59.051
     8        59.062
     9        58.896
    10        68.000
    11        59.132
    12        67.518
    13        62.715
    14        61.127
    15        60.959
    16        58.952
    17        65.498
    18        64.636
    19        64.234
    20        63.950
    21        58.740
    22        58.628
    23        58.518
    24        64.751
    25        60.742
    26        60.782
    27        60.689
    28        60.620
    29        58.782
    30        61.266
    31        61.196
    32        66.779
    33        60.347
    34        59.625
    35        67.186
    36        70.628
    37        64.413
    38        64.254
    39        61.173
    40        60.986
    41        60.859
    42        67.313
    43        66.498
    44        66.265
    45        73.208
    46        69.019
    47        59.084
    48        58.702
    49        58.575
    50        58.981
    51        58.962
    52        72.962
    53        64.662
    54        93.537
    55        74.840
    56        64.678
    57        63.313
    58        77.661
    59        81.452
    60        78.404
    61        66.582
    62        64.694
    63        59.468
    64        58.662
    65        58.550
    66        68.006
    67        80.724
    68        80.203
    69        70.724
    70        67.127
    71        58.790
    72        58.628
    73        67.547
    74        65.055
    75        63.984
    76        74.278
    77       108.038
    78        76.653
    79        74.378
    80        78.524
    81        76.770
    82        66.863
    83        65.444
    84        59.982
    85        59.075
    86        58.910
    87        67.747
    88        81.301
    89        79.237
    90        69.843
    91        69.030
    92        68.503
    93        67.855
    94        66.941
    95        77.129
    96        78.842
    97        71.306
    98        59.718
    99        69.625
   100        91.222
   101        80.421
   102        87.851
   103        68.182
   104        67.402
   105        67.289
   106        67.205
   107        67.146
   108        81.254
   109        80.770
   110        78.742
   111        66.457
   112        63.467
   113        62.836
   114        62.540
   115        62.366
   116        62.411
   117        64.203
   118        67.749
   119        66.340
   120        62.642
   121        61.011
   122        81.145
   123        77.315
   124        66.942
   125        66.077
   126        62.577
   127        62.352
   128        62.254
   129        69.026
   130        83.369
   131        80.140
   132        73.811
   133        72.626
   134        61.891
   135        68.541
   136        67.485
   137        71.349
   138        70.547
   139        73.366
   140        60.124
   141        68.149
   142        75.119
   143        81.692
   144        72.515
   145        71.287
   146        67.992
   147        86.102
   148        81.502
   149        79.439
   150        67.710
   151        66.846
   152        66.493
   153        58.712
   154        58.688
   155        61.106
   156        61.057
   157        61.079
   158        61.056
   159        61.036
   160        61.023
   161        66.400
   162        71.472
   163        59.414
   164        84.200
   165        89.496
   166        74.494
   167        63.962
   168        59.186
   169        58.593
   170        58.477
   171        58.534
   172        58.461
   173        58.419
   174        64.615
   175        60.857
   176        60.903
   177        60.767
   178        60.681
   179        72.323
   180        61.010
   181        60.896
   182        58.710
   183        67.778
   184        65.082
   185        60.712
   186        60.774
   187        65.987
   188        62.272
   189        61.394
   190        65.049
   191        64.406
   192        65.307
   193        63.143
   194        59.069
   195        58.718
   196        58.543
   197        62.174
   198        58.650
   199        58.523
   200        71.947
   201        61.087
   202        61.010
   203        66.362
   204        88.473
   205        86.344
   206        61.964
   207        61.217
   208        73.949
   209        67.857
   210        65.868
   211        69.467
   212        69.013
   213        73.596
   214        66.738
   215        59.348
   216        58.530
   217        58.399
   218        61.858
   219        58.462
   220        58.354
   221        71.973
   222        61.113
   223        61.048
   224        58.506
   225        73.831
   226        65.883
   227       106.557
   228        78.595
   229        67.672
   230        65.836
   231        69.398
   232        72.080
   233        68.893
   234        64.350
   235        64.197
   236        59.337
   237        58.518
   238        58.369
   239        64.884
   240        60.809
   241        60.890
   242        60.785
   243        60.700
   244        61.324
   245        80.348
   246        75.417
   247        80.264
   248        64.378
   249        60.345
   250        74.819

  LOWEST ENERGY:   58.354


S T E E P E S T   D E S C E N T

STEPS = 500

STEP n       E(n)         E(n-1)    
------------------------------------
    0      58.354      ----
   10    58.34134    58.34294
   20    58.33144    58.33191
   30    58.32289    58.32375
   40    58.31612    58.31645
   50    58.31017    58.31044
   60    58.30363    58.30507
   70    58.29858    58.29888
   80    58.29428    58.29417
   90    58.28891    58.29033
  100    58.28490    58.28520
  110    58.28192    58.28144
  120    58.27736    58.27826
  130    58.27401    58.27434
  140    58.27159    58.27106
  150    58.26787    58.26919
  160    58.26530    58.26557
  170    58.26402    58.26313
  180    58.26063    58.26158
  190    58.25829    58.25854
  200    58.25750    58.25657
  210    58.25463    58.25553
  220    58.25288    58.25316
  230    58.25219    58.25146
  240    58.24999    58.25099
  250    58.24866    58.24908
  260    58.24741    58.24754
  270    58.24679    58.24645
  280    58.24545    58.24623
  290    58.24427    58.24451
  300    58.24340    58.24352
  310    58.24296    58.24257
  320    58.24172    58.24218
  330    58.24095    58.24116
  340    58.24121    58.24025
  350    58.23969    58.24061
  360    58.23895    58.23933
  370    58.23824    58.23835
  380    58.23820    58.23777
  390    58.23735    58.23807
  400    58.23663    58.23683
  410    58.23616    58.23627
  420    58.23616    58.23571
  430    58.23530    58.23572
  440    58.23488    58.23506
  450    58.23546    58.23448
  460    58.23428    58.23513
  470    58.23382    58.23417
  480    58.23337    58.23347
  490    58.23359    58.23314
  500    58.23299    58.23367
1 molecule converted
23 audit log messages 

Скриптом prepare_ligand4.py из пакета Autodock tools создадим pdbqt файл лиганда:

In [1]:
%%bash
prepare_ligand4.py -l nag.pdb
Traceback (most recent call last):
  File "/usr/share/pyshared/AutoDockTools/Utilities24/prepare_ligand4.py", line 11, in <module>
    from AutoDockTools.MoleculePreparation import AD4LigandPreparation
  File "/usr/lib/python2.7/dist-packages/AutoDockTools/MoleculePreparation.py", line 11, in <module>
    import numpy.oldnumeric as Numeric, math 
ImportError: No module named oldnumeric

Скриптом prepare_receptor4.py из пакета Autodock tools создадим pdbqt файл белка, который мы получили в прошлом практикуме:

In [1]:
%%bash
prepare_receptor4.py -r ../Task9/seq.B99990001.pdb
Traceback (most recent call last):
  File "/usr/share/pyshared/AutoDockTools/Utilities24/prepare_receptor4.py", line 10, in <module>
    import MolKit.molecule
  File "/usr/lib/python2.7/dist-packages/MolKit/molecule.py", line 25, in <module>
    from mglutil.util import misc
  File "/usr/lib/python2.7/dist-packages/mglutil/util/misc.py", line 19, in <module>
    import numpy.oldnumeric as Numeric
ImportError: No module named oldnumeric

Теперь надо создать файл с параметрами докинга vina.cfg.

Для этого определим из модели seq.B99990001.pdb центр масс куба, куда будет происходить докинг, с помощью PyMol командой pseudoatom:

In [18]:
cmd.load("../Task9/seq.B99990001.pdb")
cmd.select("sele", selection="resi 120+80+62")
cmd.pseudoatom("psd", selection="sele")
cmd.save("pseudo.pdb", "psd")
In [1]:
%%bash
cat pseudo.pdb
HETATM    1 PS1  PSD P   1      40.197  41.379  28.890  0.00  0.00      PSDOPS  
END

Теперь построим файл vina.cfg:

In [21]:
%%bash
echo "center_x=40.2
center_y=41.4
center_z=28.9

size_x = 25
size_y = 25
size_z = 25

num_modes = 20 " > vina.cfg

Проведём первый докинг:

In [24]:
%%bash
vina --config vina.cfg --receptor seq.B99990001.pdbqt --ligand nag.pdbqt --out nag_prot.pdbqt --log nag_prot.log
In [33]:
cmd.reinitialize()
cmd.load("nag_prot.pdbqt")
cmd.load("../Task9/seq.B99990001.pdb")
cmd.do("set all_states, on")
cmd.show_as(representation="cartoon", selection="seq.B99990001")
Img("dok1.png")
Image(filename="dok1.png")
Out[33]:

Из файла nag_prot.log выпишем энергии 3 лучших положений и геометрическую разницу между ними:

Состояние Энергия (kcal/mol) Разница (RMSD)
1 -5.4 0.000
2 -5.2 3.408
3 -5.1 5.383

Теперь проведём докинг рассматривая подвижность некоторых боковых радикалов белка. Сначала разобьем белок на две части, подвижную и неподвижную.

Для подвижной части выберем 3 аминокислоты которые использовались в прошлом задании для позиционирования лиганда:

In [4]:
%%bash
python /usr/share/pyshared/AutoDockTools/Utilities24/prepare_flexreceptor4.py -r seq.B99990001.pdbqt -s ASN62_TRP80_ASN120
Traceback (most recent call last):
  File "/usr/share/pyshared/AutoDockTools/Utilities24/prepare_flexreceptor4.py", line 10, in <module>
    from MolKit.protein import ProteinSet, ResidueSet
  File "/usr/lib/python2.7/dist-packages/MolKit/protein.py", line 26, in <module>
    from MolKit.tree import TreeNode, TreeNodeSet, TreeNodeSetSelector
  File "/usr/lib/python2.7/dist-packages/MolKit/tree.py", line 41, in <module>
    from MolKit.listSet import ListSet
  File "/usr/lib/python2.7/dist-packages/MolKit/listSet.py", line 22, in <module>
    from mglutil.util import misc
  File "/usr/lib/python2.7/dist-packages/mglutil/util/misc.py", line 19, in <module>
    import numpy.oldnumeric as Numeric
ImportError: No module named oldnumeric

Проведём докинг:

In [4]:
%%bash
vina --config vina.cfg --receptor seq.B99990001_rigid.pdbqt --flex seq.B99990001_flex.pdbqt --ligand nag.pdbqt --out nag_prot_flex.pdbqt --log nag_prot_flex.log

Из файла nag_prot_flex.log выпишем энергии 3 лучших положений и геометрическую разницу между ними:

Состояние Энергия (kcal/mol) Разница (RMSD)
1 -5.4 0.000
2 -5.4 3.056
3 -5.1 4.138
In [17]:
cmd.reinitialize()
cmd.load("nag_prot_flex.pdbqt")
cmd.bg_color("white")
cmd.load("seq.B99990001_rigid.pdbqt")
cmd.do("set all_states, on")
#cmd.show_as(representation="cartoon", selection="seq.B99990001_rigid.pdbqt")
Img("dok2.png")
Image(filename="dok2.png")
Out[17]:

На картинке не видно, но полученное изображение отличается тем, что 3 выделенных отстатка белка приобрели подвижность.
PS к сожалению, pdbqt нельзя представить в виде cartoon

Сравним времена работы обычного и подвижного докинга:

In [21]:
%%bash
t=`time vina --config vina.cfg --receptor seq.B99990001.pdbqt --ligand nag.pdbqt --out nag_prot.pdbqt --log nag_prot.log`
echo $t

real0m7.024s
user0m44.176s
sys0m0.060s

In [22]:
%%bash
t=`time vina --config vina.cfg --receptor seq.B99990001_rigid.pdbqt --flex seq.B99990001_flex.pdbqt --ligand nag.pdbqt --out nag_prot_flex.pdbqt --log nag_prot_flex.log`
echo $t

real0m33.959s
user4m4.696s
sys0m0.156s

Время работы различается почти в 5 раз.

In []: